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SECTION  1 
INTRODUCTION 

The  modification  of  the  index  of  refraction  by  a  shock  wave  may 
produce  refraction  of  electromagnetic  signals,  thus  shock  waves  have  the 
potential  of  introducing  errors  into  radar  system  measurements  of  target 
locations.  At  high  altitudes,  i.e.,  in  the  ionosphere,  refraction  is 
determined  by  the  amount  of  ionization  present,  but  at  low  altitudes  the 
increased  electron-neutral  collisional  frequency  causes  the  accompanying 
absorption  to  dominate  any  refraction  so  produced.  Therefore  any  signifi¬ 
cant  refraction  produced  at  low  altitude  will  be  due  to  the  changes  in  the 
atmospheric  density  and  temperature.  The  atmospheric  index  of  refraction 
as  a  function  of  these  parameters  is  reviewed  in  Section  2. 

In  this  approach  to  shock  wave  induced  refraction  it  will  be 
assumed  that  the  given  parameters  include  the  actual  target,  burst  and 
emitting  source  (radar)  locations  and  that  the  resulting  angular  refrac¬ 
tive  error  is  the  desired  output.  This  is  consistent  with  the  current 
ROSCOE/NORSE  code  structure.  It  is  also  a  more  difficult  problem  to  solve 
than  predicting  the  arrival  angle  error  of  a  specified  ray,  since  the 
arrival  angle  of  the  ray  that  originates  at  the  target  is  not  a  priori 
known.  The  errors  expected  are  fortunately  relatively  small  (<  1*)  thus  a 
method  of  iteration  beginning  with  the  unrefracted  ray,  will  be  used. 
Section  3  presents  these  procedures  for  four  possible  cases,  defined  by 
the  locations  of  the  source  and  of  the  target  relative  to  the  shock 
front.  The  computer  program  and  some  preliminary  results  are  described  in 
Section  4.  Program  listings  are  included  as  an  appendix. 


Limits  on  when  refraction  due  to  spherically  expanding  shock 

waves  need  be  considered  have  been  generated.  The  early  time  limit  to  the 

importance  of  refractive  effects  is  set  by  the  cessation  of  absorptive 

effects.  Exact  values  are  of  course  sensitive  to  geometry,  however  the 

relatively  strong  dependence  of  ionization  on  temperature  limits  the 

temperature  range  of  interest  to  between  about  800*  K  for  a  long  path  to 

1200*  K  to  1400*  K  for  relatively  short  paths.  This  is  further  discussed 

in  Appendix  I  where  the  absorption  limit  is  discussed  in  detail.  An 

absolute  upper  limit  to  the  temperature  of  1500*  is  suggested.  The 

location  of  this  temperature  contour  varies  with  time.  At  the  shock 

front,  the  corresponding  shock  strength  is  about  25  psi  overpressure, 

which  occurs  at  a  scaled  range  of  56(H)1/3  meters.  The  fireball 

K  I 

however  is  significantly  hotter  than  the  shock  and  will  move  outward  to 
almost  twice  this  range. 

The  lower  limit  to  when  refraction  need  be  predicted  is  a  func¬ 
tion  of  both  the  accuracy  required  and  the  location  of  the  burst  point 
relative  to  the  sight  path.  For  example,  if  an  instantaneous  refractive 
error  of  1  milliradian  is  significant  and  if  the  shock  wave  becomes  tan¬ 
gent  to  the  sight  path  at  about  its  mid  point,  then  a  density  increase  of 
0.4  percent  is  sufficient  and  this  corresponds  to  a  shock  overpressure  of 
0.03  psi,  i.e.,  a  relatively  weak  shock. 
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SECTION  2 

INDEX  OF  REFRACTION 


The  index  of  refraction  at  any  point  is  a  function  of  both  the 
amount  of  ionization  present  and  of  the  density  and  temperature  of  the 
various  neutral  constituents.  In  analyzing  refractive  errors  for  propaga¬ 
tion  paths  at  low  altitude,  it  is  only  necessary  to  consider  the  changes 
to  the  neutral  constituents.  This  is  because  the  absorption  of 
electromagnetic  energy  is  also  proportional  to  the  amount  of  ionization 
present  and  to  the  electron-neutral  collisional  frequency,  which  is 
proportional  to  pressure.  Therefore  at  the  lower  altitudes  high  levels  of 
absorption  will  occur  on  typical  propagation  paths  that  have  even  moderate 
amounts  of  refraction.  This  can  be  seen  by  comparing  the  angular 
deviation  of  the  ray  path,  in  degrees,  to  the  one-way  absorption,  A  in 
dB.  Reference  1  gives  the  approximate  relation*  for  the  case  of 
spherically  stratified  ionization 


j.  _  107 
A  v  r 


7db 


(2-1) 


where  v  is  the  electron-neutral  collision  frequency  (sec-1)  in  the  region 
and  r  is  a  characteristic  dimension  in  (km).  Using  typical  values,  v  « 
1011  sec*1  and  r  «•  ~  km,  we  obtain  ^  -  10*3  degree  per  db  or  16  micro¬ 
radians  per  db,  i.e.,  the  absorption  along  a  path  that  yields  1  milli- 
radian  deflection  due  to  an  ionization  gradient  will  also  yield  60  db  of 


See  for  example  Page  6-6  of  Reference  1. 
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absorption.  In  the  Appendix  it  is  shown  that  when  the  temperature  exceeds 
about  1000*  K  -  1500*  K,  then  the  quasi-equilibriun  ionization  resulting 
from  delayed  fission-product  radiations  will  produce  such  high  levels  of 
absorption.  Therefore  when  predicting  refraction  caused  by  a  shock  wave 
it  is  only  necessary  to  consider  the  impact  of  density  and  temperature 
changes  on  the  index  of  refraction,  and  only  for  temperatures  below 
1500*K. 


A  best  fit  to  the  data  on  the  atmospheric  refractive  index  at 
radio  frequencies  was  determined  at  NBS  and  reported  in  Reference  2  to  be 
(their  Equation  7) 


N  =  (n-l)106  =  77.6  j  -  6  i  +  3.75  x  10s  (2-2) 

where  P  is  the  total  pressure,  in  millibars,  T  is  the  absolute  tempera¬ 
ture,  in  *K,  and  e  is  the  partial  pressure  of  water  vapor,  in  millibars. 

A  form  that  contains  the  wavelength  dependence  can  be  used  to  show  that 
this  dependence  is  negligible  at  radio  frequencies.  This  form  is  as  given 
by  Allen  (Reference  3)  at  STP  as 

(n-1)  x  106  -  64.328  +  — 2?4.9_  l.  +  2-55,4  ,  at  15’C  (2-3) 

146  -  (|)2  41  -  (j)2 

where  x  is  the  vacuum  wavelength  in  microns.  Since  we  are  not  interested 
in  wavelengths  less  than  one  mm  or  103  microns,  the  corrections  are  always 
negligible  and  this  reduces  to  (n-1)  x  106  =  272.6,  consistent  with 
Equation  2-2  above. 

It  is  convenient  to  simplify  Equation  2-2.  We  introduce  the 
relative  partial  pressure  of  water  vapor,  e'  =  c/p.  The  gas  law  allows 


p/T  to  be  replaced  by  the  density,  p,  in  cpi/cm3,  which  then  can  be 
factored  out.  Following  Reference  3,  we  can  also  combine  the  second  term 
with  the  third  term  with  only  small  error  for  the  temperature  range  of 
interest  -  especially  since  the  ambient  value  of  e'  is  highly  variable. 


N  =  (n-l)106  =  2.2  x  105p  (1  +  4.8  x  103  £-) 


<2-4a) 


n  =  1  +  0.22p  (1  +  4.8  x  103  1_)  . 

T 


(2-4b) 


Later  it  will  be  convenient  to  consider  the  ratio  of  the  indices 
of  refraction  across  a  shock  front.  Using  the  subscripts  o  and  s  to 
signify  the  ambient  and  shocked  conditions,  this  ratio  is 


n  1  +  0.22  p  (1  +  4.8  x  103  *  ) 

_i  =  s  ^ _ S.  (2-5) 

n0  1  ♦  0.22  p0  (1  t  4.8  x  103  £l) 


which  can  be  closely  approximated  by 

P  0  i  0  * 

-1  -  [1  +  0.22  p$  (1  +  4.8  x  103  _)]  [1  -  0.22  Pq  (1  +  4.8  x  103  — )] 


1  +  0.22  p  [(II  -  1)  +  4.8  x  103  11  [J.  •  _£  -  1)]  . 
0  o  T  p^  T 


(2-6) 


o  po 


In  this  form  we  see  that  an  accurate  knowledge  of  the  water  vapor  content 
is  usually  not  important.  In  much  of  the  range  of  interest  the  compres¬ 
sion  ratio  is  about  40%  greater  than  the  shock  temperature  ratio. 
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Using  an  ambient  temperature  of  288*  K  the  second  term  in  the  bracket 
becomes  6.7  e*.  Typical  values  of  e1  are  less  than  0.03,  yielding  a  vali.a 
of  0.2  or  less,  to  be  compared  with  the  shock  overdensity  ratio,  (ps/p0-l) 
=  (p-1),  which  is  generally  much  greater. 

For  discussion  purposes,  it  is  convenient  to  use  an  ambient 
density  of  10-3  g/cm3  (and  a  value  corresponding  to  an  altitude  of  about 
6000  ft  above  sea  level)  and  a  relatively  high  value  of  water  vapor 
partial  pressure  of  3%  (corresponding  to  saturated  air  at  about  25*  C). 
Equation  2.6  then  becomes 

nc  1  Tn 

—  =  1  +  2.2  x  10--  [(u-1)  *  -  (m  -  1)]  .  (2-7) 

*0  2  Ts 

It  is  shown  in  the  Appendix  that  if  the  shock  temperature  exceeds  about  4 

times  ambient  temperature  then  absorption  is  the  dominant  process.  The 

corresponding  compression  ratio,  u  =  £i,  is  4.7.  Inserting  these  values 

PO 

into  Equation  2-7  yields  an  upper  limiting  value  of 

—  =  1  +  2.2  x  10-**  [3.7  +  I  (ll l  -  1)]  =  1  +  8.3  x  10-4  (2-8) 

n„  2  4 
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SECTION  3 

SHOCK-INDUCED  REFRACTIVE  ERRORS 


The  geometry  of  the  shock  wave  related  refraction  problem  is 
shown  in  Figures  1  and  2.  The  ambient  (external)  atmosphere  is  assumed  to 
be  uniform  and  the  shock  wave  is  assumed  to  be  spherical.  Therefore  the 
source,  target  and  burst  points  define  the  plane  of  these  figures.  There 
is  no  refraction  out  of  this  plane.  Figure  1  shows  the  larger  view  of 
this  plane.  A  source  (e.g.,  radar)  is  located  at  the  bottom  of  the 
figure,  tracking  a  target  at  the  top  which  is  moving  to  the  left.  A  burst 
occurs  at  a  range  Rg  from  the  source.  The  angle  at  the  source  between  the 
sight  lines  to  the  burst  and  to  the  real  location  of  the  target  is  6}.  A 
shock  wave  is  expanding  from  the  burst,  at  a  radius  Sp.  The  density 
increase  within  the  shock  wave  produces  an  increase  in  the  local  index  of 
refraction.  Refractive  effects  will  occur  only  after  the  target  passes 
behind  the  shock  front,  i.e.,  not  until  after  it  passes  point  T0.  When 
viewed  through  the  shocked  region  the  apparent  target  location  T1  will  be 
to  the  right  (in  this  figure)  of  the  real  location,  T.  As  the  shock  wave 
intercepts  the  line  of  sight  there  will  be  a  portion  of  the  trajectory 
that  is  not  visible.  Of  course,  in  the  general  case  the  current  shock 
location  could  be  beyond  either  the  target  or  source  location  or  beyond 
both. 


Figure  2  illustrates  a  possible  propagation  path  within  the 
shocked  region.  The  insert  at  the  top  of  the  figure  shows  a  possible 
radial  profile  of  density.  This  profile,  including  the  current  shock 
radius  is  assumed  to  be  given.  For  example  it  may  be  obtained  from  the 
Nuclear  Blast  Standard  (1  KT)  (Reference  4),  the  LAMB  code  (which  includes 
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Figure  1.  Geometric  definitions  outside  the  shock. 
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the  1  KT  standard)  or  other  simple  fits  to  the  well  known  blast  wave  data. 
In  Figure  2  the  reference  axis  has  been  chosen  along  the  radial  between 
the  burst  point  and  the  entry  point  of  the  ray  path.  Note  that  the  loca¬ 
tion  of  this  entry  point,  and  thus  also  the  angle  <p0  are  to  be  determined, 
and  when  known  represents  the  solution  of  the  problem. 

The  location  of  the  target,  labeled  T,  may  be  anywhere  along  the 
path,  including  inside  the  shock.  The  apparent  location,  T',  will  be 
along  the  dashed  extension  of  the  incident  ray,  at  different  range.  The 
propagation  path  starts  from  S  at  the  angle  (e1  +  59J  and  enters  the 
shock  at  point  A  where  it  makes  the  angle  ^  with  a  radial  from  the  burst 
point.  The  actual  ray  to  the  target  is  refracted  towards  this  radial 
through  an  angle  Sip.  When  both  the  radar  and  target  are  outside  the 
shocked  region  the  ray  path  exits  the  shock  at  point  C  where  it  is  refrac¬ 
ted  away  from  the  radial  through  C  by  the  same  angle  6ip,  and  arrives  at 
the  target  location  at  an  angle  692  off  of  the  line  S-T.  In  this  case  of 
spherical  symmetry  it  is  possible  to  replace  the  actual  curved  path  by  a 
straight  line  connecting  the  entry  and  exit  points.  This  is  equivalent  to 
approximating  the  radially  varying  value  of  the  index  of  refraction  by 
some  'effective1  value,  n,  yet  to  be  defined. 

REFRACTION  IN  A  SPHERICALLY  STRATIFIED  REGION 

Referring  to  Figure  2,  the  polar  coordinates  of  an  arbitrary 
point  along  the  ray  path  interior  to  the  shock  as  measured  from  the  burst 
point  are  (r,  <p),  and  as  measured  in  polar  coordinates  from  the  entry 
point  are  (p,  ip).  The  angle  between  the  vector  r  and  the  local  direction 
of  the  ray  path  is  q.  Note  from  the  geometry  of  these  definitions  that 

p  sin  ip  •=  -r  sin  <p 

and 

p  cos  ip  -  SR-r  cos  <j> 
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thus 


y  -  tan"1 


r  sin  4> 
'Sp-r  cos  <fi 


)  • 


(3-1) 


The  change  in  ray  direction  upon  shock  entry  (and  when  appropri¬ 
ate  upon  exit)  is  obtained  from  Snell's  law.  Note  that  at  entry  point 
Thus 


n0  sin  *0  =  sin  ^  =  n$  sin  .  (3-2) 

The  solution  for  the  ray  path  within  a  spherically  stratified  region  has 
been  given  by  Archer  (References  5  and  6).  The  form  of  Snell's  law  for  a 
spherically  stratified  medium  is 

r  n(r)  sin  c  =  K  (3-3) 

The  radial  variation  of  n  is  obtained  by  substituting  the  given  radial 

variations  of  density,  p(r),  and  of  temperature,  T(r),  into  equation 

(2-4b)  of  the  previous  section.  The  constant  K  may  be  defined  by  values 

either  upon  entry  or  at  the  point  of  closest  approach  where,  c  =  n/2  and  r 

=  R  .  Combining  equations  3-2  and  3-3  provides  a  relation  between  R 
m  m 

and  the  (as  yet  unknown)  angle  of  incidence  outside  the  shock,  yields 


=  r  n(r)  sin  x. 

( 3-4  a ) 

-  Rm 

(3-4b) 

=  SR  n$sin  x,%  *  SR  sin  ^ 

(3-4c) 

•  SR  nQ  sin 

( 3— 4d ) 
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The  differential  equation  of  the  ray  path  within  the  shocked 


region  is 


+  K  dr 


/r2  n2(r)  -  K2 


(3-5) 


The  integral  of  this  equation  between  two  points  (r  1  ,<^>1 )  and  (r2,  412) 
provides  the  angle  $  =  $2  -  $1  between  these  points,  i.e.. 


4> 


_ K_dr _ 

r  /r2  n2 (r)-K2 


(3-6) 


where  we  have  used  the  symbol  (p  to  indicate  the  integral  is  to  be  taken 
along  a  specific  path  and  not  just  between  the  radial  points  rx  and  r2. 
Specif ical ly,  for  paths  which  pass  through  the  point  of  closest  approach 
this  integral  runs  from  rl  to  Rm  to  r2.  In  the  following,  the  point  rt  is 
either  the  known  location  of  the  source  when  the  source  is  within  the 
shock  or  the  entry  point  of  the  ray  path  when  the  source  is  outside.  Note 
that  while  the  radius  to  that  entry  point  is  known  (i.e.,  SD)  its  angular 
position  is  not  known.  Similarly  r2  refers  to  either  the  target  location 
or  the  exit  point  of  the  ray  path. 


Unfortunately  the  constant  K  in  equation  6  contains  Rm  which 
depends  on  the  ray  path,  or  equivalently  the  angle  i|>0  which  depends  on  the 
entry  point.  Therefore  it  is  necessary  to  solve  the  set  of  equations 
iteratively.  Fortunately,  we  are  interested  only  when  the  refractive 
error  is  small,  thus  the  undeviated  path  provides  a  first  estimate  of  R^ 
or  i|>0.  From  Figure  1  this  can  be  seen  to  be 

Rm  “  RB  Sln  81  (3'7) 
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Using  this  value  of  r,  a  density  and  temperature  are  obtained  from  the 
given  shock  wave  profile.  These  are  then  used  in  equation  2- 4b  to  obtain 
the  local  index  of  refraction.  This  together  with  the  radius  then 
provides  the  initial  estimate  of  K. 

BEARING  ERROR  PREDICTIONS 

There  are  four  possible  geometric  situations  that  may  occur 

depending  on  the  relative  location  of  the  shock  to  the  source  and  to  the 

target.  Each  requires  a  slightly  different  methoo  of  solution,  although 

the  principle  remains  the  same.  Using  the  assumption  that  the  refractive 

error  is  small,  K  as  obtained  above  is  used  to  determine  a  first  estimate 

of  the  angular  extent  of  the  ray  path  interior  to  the  shock,  4.  In  three 

of  these  cases,  geometric  considerations  then  provide  an  estimate  of  the 

shock  front  entry  or  exit  angle,  which  provides  a  second  estimate  of  K  and 

leads  to  a  situation  that  may  be  iterated.  In  the  fourth  case,  which  we 

shall  discuss  first,  the  iteration  is  on  the  location  of  R  since  the 

m 

shock  front  is  never  reached. 

Source  and  Target  Both  Within  the  Shocked  Zone 

In  this  case  the  angular  extent,  $,  must  be  equal  to  the  angle 
at  the  burst  point  between  the  radii  to  the  source  and  the  target,  i.e., 

♦  *  «  -  (0i+e2)  (3-8) 

The  end  points  of  the  integral  are  also  determined  by  geometry,  thus  the 
only  available  parameter  is  K,  or  equivalently  the  location  of  the  minimum 
radius  to  the  propagation  path.  Given  that  parameter  the  initial  direcion 
of  the  ray  is  obtained  from  equation  3-4a.  The  refractive  error  is  then 
the  difference  between  this  angle  and  0i,  i.e., 
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(3-9a) 


«ej  *  C(Rb)  -  0! 


«  sin 


.1 


Rb  n(R g) 


(3- 9b) 


An  estimate  of  the  correction  to  K  necessary  to  make  the  calcu¬ 
lated  value  of  $  agree  with  the  geometric  value  may  be  obtained  by  differ¬ 
entiating  equation  3-6  w.r.t  k. 


d£  _  d 
dk  dk 


0  —=z - +  F(r2,K)  ^  -  F(rlfK)  i£l 

rl  r  /rV(r)  -  K2  K 


(1  + 


+  0-0 


r  n2(r)  dr 
(r2n2(r)  -  K2)3/2 


(3-10) 


Source  Outside  and  Target  Inside 

This  situation  is  illustrated  in  Figure  3.  The  initial  estimate 

of  6  is  based  on  K  as  obtained  from  R  ,  equation  3-4b.  This  then  defines 

m 

a  revised  entry  point  and  thus  a  new  value  for  from  geometry,  i.e.,  two 
sidfes  (R B  and  SR)  and  the  included  angle  ( 93-4>)  of  the  triangle  S-B-0  are 
known.  The  revised  value  of  K,  obtained  from  equation  3-4d,  is  then  used 
to  recalculate  <(>•  The  process  is  iterated  for  a  consistent  set  U,\|)0)- 
The  refractive  error,  by  geometry  is  then 


60!  =  ♦o"0i”(e3"+)  *  'ko+02"1,+<1 


(3-  ’) 
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Figure  3.  Source  outside  and  target  inside  the  shocked  region. 

Note  that  the  curved  nature  of  the  path  interior  to  the  shock  does  not 
enter  directly. 

Source  Inside,  Target  Outside 

This  is  essentially  the  inverse  of  the  previous  tase  with  one 
major  difference,  since  the  refractive  error  occurs  on  the  curved  portion 
of  the  path  which  extends  to  the  source.  The  method  of  solution  is  essen¬ 
tially  the  same,  except  that  once  the  consistent  set  U,\|>0)  is  obtained 
the  corresponding  value  of  K  is  used  in  equation  3-4a  to  obtain  the 
instantaneous  ray  path  angle  5  at  the  radius  RD.  The  refractive 

O 

error  is  then 

501  ~  5  ( Rg )  “  ®i 

which  is  the  same  as  equation  3-9a. 
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Target  and  Source  Both  Outside 

When  the  propagation  path  completely  traverses  the  shocked 
region  as  depicted  in  Figures  1  and  2,  then  the  radial  from  B  to  D  will  be 
perpendicular  to  this  path  and  the  angle  4»0  (or  B-A-C)  is  the  compliment 
of  the  angle  <k  (or  A-B-D),  where  $  is  defined  as  one-half  the  total 
angular  extent  of  the  ray  path  inside  the  shock,  i.e., 

♦m  '  f"  — --K— .  (3-12) 

SR  r  /rV(r)  -  K2 

The  effective  value  of  the  index  of  refraction  which  corresponds  to  a 
straight  line  path  inside  the  shock,  is 

ne  -  n0  sin  4>0/sin  4^  =  n0  sin  4»0/cos  (3-13) 

Consider  the  triangles  formed  by  the  unrefracted  extension  of  the  external 
ray  paths.  From  the  small  triangle  the  angle  between  lines  S-A-D  and 
D-C-T  is  264;.  From  t,ie  large  triangle  (S,T,D)  that  angle  is  60!  +  602. 
Thus 

50j  +  602  =  2(<l>o"’l'e)  =  2641  (3-14) 

The  angle  between  the  radial  BO  and  the  radial  that  is  perpendicular  to 
ST  can  be  seen  to  be  (54,-50!)  which  is  also  1  ( 602- 60 1 ) . 

Applying  the  law  of  sines  to  the  triangles  S-E-D  and  D-E-T 

yields 


(Ri+d)  sin  =  t  sin(90-«*)  =  (R2-d)  sin  502  ,  (3-15) 
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where 


Ri  =  Rg  cos  0! 
R2  *  Ry  cos  02 


60o  _tfl 

d  =  Rg  sin  6i  •  tan( - — -i) 

691  and  692  are  small  angles,  thus  the  sines  in  equation  3-15  may  be 
approximated  by  the  angles,  yielding 

( R ! +d )  69 L  -  (R2-d)  fie2  . 

Using  this  in  Equation  3-14  to  eliminate  se2  yields 
U  +  501  *  264-  . 

R2  “0 

d  can  also  be  neglected,  being  much  smaller  than  either  Ri  or  R2,  yielding 
(1  +  J?i)  sqj  =  2Si/>  .  (3-16) 


A  relation  between  the  viewing  angle,  0 1 +60 j  and  the  angle  of  incidence, 
of  the  ray  SA  at  the  shock  front  is  obtained  by  applying  the  law  of 
sines  to  triangle  S-B-A. 


We  now  apply  Snell's  law  in  the  form 


n0  sin  4»0  =  ne  sin  (*0-6i|>)  (3-18) 

where  ne  is  the  "effective"  value  behind  the  shock.  Expanding  the  left 
side  of  Equations  3-18  and  rearranging  yields 

,  sin  6ij> 

tan  *o  =  _ _  (3-19) 

cos  Sij>  -  no/ne 


To  eliminate  \p0  we  rewrite  sin  4>0  in  Equation  3-17  in  terms  of  the  tangent 
.  and  substitute  from  Equation  3-19,  thus 


R 


B 


—  sin  (ej+seJ  =  sin  i|>0  « 
SR 


tan  fo 


sin  Sip 

/sin2  +  cos25^  -  2  (n0/ne)cos  6i|>  +  (n0/ne)2 


sin  6ip 

+  2  ZlQ.  (l  -  cos  6<|>) 


(3-20) 


Equation  3-16  is  used  to  eliminate  5^  yielding  an  equation  containing  60! 
and  known  quantities. 


"r 

—  sin  (e j^+50 x )  - 
SR 


Sln  59^ 


'(1  -  — )2  +  2  [1  -  COS  (*1*2)  691] 


(3-21) 


2R2 


20 


We  expect  60!  to  be  a  small  angle  and  do  not  expect  the  factor  )  to 

2R  2 

prevent  a  small  angle  approximation,  then 


Rb 

_  (sin  0!  +  60i  cos  0 i ) 
SR 


(RCR,2.)  601 

iR2 - -  (3-22) 

/( i  ™  ulEE) 

n  n  2R2 
e  e 


From  Equation  2-8  we  expect  the  first  term  in  the  square  root  to  be  of  the 
order  of,  but  less  than  (8.3  x  10'1*)2.  This  generally  will  be  smaller 
than  the  second  term  whenever  refraction  exceeds  a  few  milliradians, 
therefore  as  a  first  estimate  we  neglect  the  first  term,  yielding 


or 


Rr 

(sin  +  60J 
SR 


COS  01  )  - 


601 


COS  @1 


(3-23) 


(3-24) 


REAPPEARANCE 


As  the  line  of  sight  to  the  target  passes  behind  the  shock  front 
the  target  will  disappear  (ignoring  diffraction)  and  will  not  reappear 
until  the  external  rays  S-A  and  C-T2  are  tangent  to  the  shock,  i.e.,  the 
angle  of  incidence,  ♦g  is  90*.  From  equation  3-18  and  3-16, 


.HEL  =  sin  (-  -  Sij;)  -  cos  6<ji 


*  cos 


R1+R2 

2P2 
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(3-25) 


Again  in  the  small  angle  approximation 


SSL «  l  - 1  se,]2  (3-26) 

ne  2  1  2R2  1J 

Inverting, 

60i  =  -i-1-  *4(1-  ill)  (3-27) 

Ri+R2  n$ 


Substituting  the  index  of  refraction  ratio  from  equation  (2-6)  gives 


«el  =  M.22)  p0[(u-1)  +  4.8  x  103ll  (plfl-  1)]  ,  (3- 

Rl+K2  t0  t 


28) 


for  p0  »  10" 3 


op  ,  _ i  T  1/2 

601  =  ( .02)  [(u-1)  +  4.8  x  103  1-  (v  III  -  1)  ]  ,  (3-29) 

Ki+W2  t0  i$ 

when  the  shock  first  becomes  transparent  at  Tg  ■  4  T0  and  u  *  4.7,  the 
square  root  term  yields  about  a  factor  of  2  thus  the  refraction  error  will 
be  about  40  mill iradians  times  the  factor  involving  the  relative  locations 
(which  is  usually  of  the  order  of  unity). 

If  we  define  the  minimum  value  of  interest  for  69j  we  can  then 

determine  the  minimum  shock  strength  andtherefore  the  maximum  shock  radius 

of  interest.  We  can  convert  equation  3-29  into  terms  of  the  relative 

shock  overpressure,  it  s  —  -  1,  using  the  Hugoniot  relation 

Po 


7  *  6  it 
7  +  it 
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and  the  ideal  gas  law 


lo  =  £i  .  £_  =  JL 

Ts  p  p0  1+ir 


Substitution  of  these  into  equation  (3-29)  yields 


60!  »  (-— 2--1  (.02)  \—  +  JLil  (21  *  -  ^L)]1/2 

Rl+R2‘  7+ir  (7+tT)T  1  1+w  J 


which  for  weak  shocks  reduces  to 


60i 


(_2_ 2  .)  (.02)  (1+lOe1)]1/2  ,  *  «  1 

Rj+R2  7 


■  5it 


(3-30) 


(3-31) 


For  example,  to  produce  the  above  estimate  to,  say,  1  milliradian,  would 

require  that  the  factor  within  the  square  root  be  only  — which  occurs  at 

400 

about  an  overpressure  of  about  0.03  psi,  i.e.,  in  the  far  field  of  the 
shock  wave. 


SECTION  4 

A  MODULE  TO  PREDICT  REFRACTION 
AND  PRELIMINARY  RESULTS 


A  set. of  subroutines  was  prepared  with  the  intent  that  an  appro¬ 
priate  subset  could  easily  be  adopted  into  larger  programs.  These  sub¬ 
routines  and  their  functions  are  described  in  the  following  paragraphs. 

Program  BLAST  and  subroutine  GETINPUT  provide  our  stand-alone 
driver  that  would  be  replaced  by  a  calling  procedure  within  the  larger 
program.  Inputs  that  are  expected,  units  used,  and  where  appropriate 
default  values,  are: 

1)  RHOA  =  Ambient  air  density  {gm/cm3,  default  -  1 . 225 x 10~ 3 ) 

2)  TEMPA  =  Ambient  air  temperature  ("K,  default  £  288) 

3)  WATER  =  Relative  partial  water  vapor  pressure  (default  =  .01) 

4)  W  =  Effective  blast  yield  (kilotons) 

5)  TIME  =  Time  of  interest  after  burst  (sec) 

6)  RBS  -  Actual  range  between  source  and  burst  (cm) 

7)  RTS  =  Actual  range  between  source  and  target  (cm) 

8)  THETA1  =  Actual  angle  between  RBS  and  RTS  (radians) 

In  the  stand-alone  version  these  are  obtained  via  a  common  block  from 
subroutine  GETINPUT  in  which  the  default  values  are  stored. 

SUBROUTINE  REFRACT 

Subroutine  REFRACT  is  the  heart  of  the  calculational  procedure. 
The  first  step  is  to  scale  ranges  and  time  to  equivalent  one  kiloton 


values  so  that  the  shock  location  (i.e.,  radius)  and  density  profile  data 
can  be  obtained  form  the  AFWL  Nuclear  Blast  Standard  (1  KT)  (Reference  4) 
via  a  call  to  DENSITY.  Then  a  set  of  tests  is  performed  to  determine  the 
relative  locations  of  the  source  and  the  target  relative  to  the  shock 
front  and  whether  or  not  the  LOS  extends  to  the  point  of  closest  approach 
(RMIN)  of  the  LOS  to  the  burst.  These  tests  define  the  approximation 
procedures  and  the  integration  limits. 

The  primary  output  provided  is  the  refractive  error,  in  radians, 
in  the  plane  formed  by  the  burst  point,  target  and  source:  that  is,  the 
increase  in  the  angle  THETA1  caused  by  the  refraction  of  the  ray  path  in 
passing  through  the  density  profile  of  the  blast  wave.  Additional  outputs 
that  are  available  include  PSIO,  the  angle  of  incident  at  the  shock  front 
(when  the  source  is  outside);  relative  locations  of  target  or  source  with 
respect  to  the  shock,  and  the  interior  angle  4;  and  n,  the  effective  index 
of  refraction  (when  both  target  and  source  are  outside). 

SUBROUTINE  ETA 

Subroutine  ETA  returns  the  index  of  refraction  within  the 
shocked  region  according  to  Equation  2-4b.  The  input  is  the  scaled  radial 
dimension  of  the  point  of  interest.  This  subroutine  then  calls  the  1  KT 
blast  model  to  obtain  the  overdensity  via  the  common  block  /WFRT/.  A 
temperature  is  needed  in  Equation  2-4b  as  part  of  the  water  vapor  correc¬ 
tion.  The  1  KT  blast  model  does  not  provide  a  temperature  or  internal 
energy  profile  behind  the  blast  wave.  Although  a  temperature  could  be 
obtained  from  the  overpressure  profile  and  the  equation  of  state  through 
an  iteration  procedure,  a  simple  approximation  has  been  used  instead.  The 
temperature  is  estimated  by  assuming  a  gamma  law  expansion  (at  y  -  1.4) 
from  the  current  shock  front  density,  but  limited  to  be  at  least  ambient 
temperature.  This  is  rationalized  as  being  sufficient  since  the  shock 
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temperature  does  not  exceed  1200*  K,  thus  any  error  in  the  temperature 
will  probably  be  less  than  the  uncertainty  in  the  partial  pressure  of  the 
water  vapor. 

SUBROUTINE  INTEGRT 


This  subroutine  performs  the  integration  of  Equation  3-6  to 
provide  the  angle  <j>  between  radial  limits  supplied  to  it  by  subroutine 
REFRACT.  The  procedure  used  is  to  subdivide  this  interval  into  steps 
within  which  the  value  of  eta  is  essentially  constant.  Equation  3-6  can 
then  be  integrated  analytically  to  give  the  increment,  ,  over  each  step, 
Ar,  i.e.. 


where  n.  is  the  value  of  the  index  of  refraction  at  the  center  of  the  ith 
interval.  The  approximate  number  of  steps  to  be  taken  is  specified  by  a 
data  statement  as  NUM,  however  the  step  size  can  be  decreased  or  increased 
internally  based  on  a  test  of  the  relative  change  of  the  index  of 
refraction  within  each  step. 

The  calling  sequence  for  this  subroutine  expects 


RS  =  First  integration  limit  of  scaled  radius 

RM  =  Second  integration  limit 

CK  *  Constant  K  of  equation  3-6 

ETANEW  =  Value  of  eta  at  starting  point  of  integration 

DPDK  =  A  trigger  value  which  causes  a  calculation  of  ~  when 

dK 

positive  or  zero,  but  skips  this  calculation  when  set 
negative 


Upon  return,  in  addition  to  providing  the  values  of  4*  and  when 

requested  the  subroutine  stores  the  most  recent  values  of  n.  in  ETANEW 
dK  1 

for  use  in  the  second  integration  when  needed. 

SUBROUTINES  DENSITY  AND  AIRPT 

Subroutine  DENSITY  contains  only  those  portions  of  the  1  KT 
blast  standard  that  are  required  for  the  refraction  prediction.  These 
were  extracted  from  Reference  4  and  are  carried  as  a  separate  routine  so 
that  when  the  more  complete  set  of  blast  subroutines  is  used  elsewhere 
within  a  larger  code,  this  subroutine  can  be  deleted  and  that  set  used. 

An  initial  call  to  DENSITY  (TIME)  at  each  time  of  interest  sets  the 
following  parameters. 

PRAD  =  Shock  front  radius  (cm) 

OPPK  -  Peak  overpressure,  at  PRAD  (dyne/cm2) 

ODPK  =  Peak  overdensity,  at  PRAD  (gm/cm3) 

RDZ  =  Radius  at  which  overdensity  passes  through  zero  (cm) 
TEMPK  =  Shock  front  temperature  (*K) 

Subsequent  calls  at  the  same  time  use  entry  DENS(RAD)  and  obtain  ODR,  the 
overdensity  at  the  specified  radius,  RAD.  Outputs  from  DENSITY  are  trans¬ 
ferred  via  the  common  block  /WFRT/.  Subroutine  DENSITY  requires  the  air 
equation  of  state  to  calculate  ODPK,  from  the  prediction  value  of  OPPK. 

The  DOAN-NICKLE  equation  of  state  of  air  as  given  in  subroutine 
AIRPT(E,R,G,P,T)  of  the  MDAC  version  of  LAMB  has  been  used,  since  it 
includes  the  temperature  and  pressure  thus  providing  TEMPK.  This  subrou¬ 
tine  appears  to  include  the  subroutine  AIR(E,R,G)  of  the  1  KT  standard, 
which  may  be  accessed  through  an  entry  call. 
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TYPICAL  RESULTS 


These  subroutines  have  been  exercised  for  the  several  conceiv¬ 
able  types  of  sight  paths,  depending  on  the  relative  location  sof  the 
radar  target,  shock  wave  and  point  of  closest  approach  to  the  burst  point 
of  the  sight  path.  For  those  case  in  which  absorption  dominates  or  where 
the  shock  is  not  intersected  a  message  to  that  effect  is  produced  without 
the  prediction  of  refraction.  Table  1  describes  eight  possible  sight 
paths  through  the  shock  wave  from  a  one  kiloton  burst  at  one  second  and 
lists  the  calculated  refractive  error. 

Table  1.  Sight  path  parameters  and  results. 


Case 

RBS 

RTS 

01 

60! 

cm 

cm 

radians 

mil  1  iradi ans 

1 

1.3E5 

IE5 

.1 

8.4E-4 

2 

1.3E5 

1.6E5 

.15 

2.9E-2 

3 

1.3E4 

2E4 

.7 

2.5E-4 

4 

3E4 

5E4 

.72 

-6.0E-5 

5 

2.5E4 

3E4 

1.9 

658. (?) 

6 

4E4 

2E5 

.7 

0.136 

7 

3E4 

1.5E5 

2.1 

-1058. (?) 

8 

1.5E5 

2.5E5 

.22 

16  to  104(*) 

(?)  These  values  are  abnormally  large;  indicating  that  the  "small 

derivation"  approximation  is  invalid  but  do  show  that  refraction  will 
be  an  extremely  severe  problem. 

(*)  This  case  did  not  converge  but  oscillated  between  these  values,  again 
indicating  a  serious  refractive  error. 
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APPENDIX  I 
ABSORPTION 


When  the  electron  concentration  within  a  shock  wave  is  suffi¬ 
ciently  high,  the  absorption  of  electromagnetic  waves  becomes  so  great 
that  refraction  effects  can  be  ignored.  Here  we  will  generate  estimates 
of  the  conditions  at  which  such  circumstances  occur  and  identify  those 
cases  wherein  refraction  might  be  important.  Except  where  noted,  the 
basic  equation  used  below  are  from  The  Aids  for  the  Study  of  Electromag¬ 
netic  Blackout  (Reference  1). 

The  differential  absorption  of  an  electromagnetic  wave  of 
angular  frequency,  u  (radians/sec) ,  can  le  expressed  as 


Ab 


46  gv 

(gv)2  +  (hu)2 


db/m  (1-1) 


where  Ng  is  the  local  electron  concentration  (e/cm3)  and  v  is  the  electron 
collision  frequency  (sec'l).  At  altitudes  below  about  100  km  the  colli¬ 
sion  frequency  of  importance  is  that  with  neutral  particles,  which  is 

v  -  1.7  x  105  p  ,  sec-1  (1-2) 

where  p  is  the  local  air  pressure  (dynes/cm2).  At  sea  level  the  ambient 
pressure  is  pQ  *  106  dyne/cm2  hence  vQ  -  1.7xl011  sec-1  which  corresponds 

to  an  operating  frequency,  f  =  —  =  —  =  2.7xl010  -  27  Ghz.  The  operating 

2v  2n 

frequencies  of  interest  may  be  on  either  side  of  the  collision  frequency. 
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The  factors  g  and  h  in  Equation  1-1  above  are  correction  factors 
to  account  for  the  velocity  dependence  of  electron-neutral  collisions. 

For  representative  calculations  values  of  these  factors,  as  obtained  from 
graphs  in  Reference  1,  are  given  in  Table  1-1  below 


Table  1-1.  Typical  values  of  g  and  h. 


g 

h 

It)  <  V 

0.6  to  0.65 

>  1.7 

it)  **  .4  v 

0.75 

1.3 

It)  >  V 

1 

1 

For  on  <  v  absorption  is  independent  of  the  operating  frequency 
and  Equation  (1-1)  can  be  written 

Ab  -  75  ,  u  €  v,  db/m  (1-3) 

v 

Using  the  sea  level,  ambient  value  of  v  this  becomes 

Ab  =  4.5  X  10-10  Ne,  f  «  27  Ghz,  db/m  (I-3a) 

At  10  Ghz  Equation  (1-1)  becomes 

Ab  =  61  Ne _  ,  f  =  10  Ghz,  db/m  (1-4) 

v  1  +  (.64  El) 

P  . 
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Using  sea-level  conditions  in  Equation  1-4  yields 


Ab  -  2.5  x  10- 10  Ne  ,  db/m  .  (I-4a) 

Note  that  the  second  factor  in  the  denominator  of  Equation  1-4  reduced  the 
absorption  by  about  40%,  but  will  contribute  much  less  in  a  strongly 
shocked  region. 

For  frequencies  well  above  the  collision  frequency,  i.e,  milli¬ 
meter  waves.  Equation  (1-1)  becomes 

Ab  =  46  —  N  ,  u,  »  v  db/m  ( 1-5) 

a.2  e 

Probably  the  minimum  path  length  through  the  shocked  region  that 
is  ever  of  interest  will  be  a  few  tens  of  meters.  Normally  the  path  will 
be  much  longer,  i.e.,  one-hundred  or  more  meters.  Equation  I -3a  indicates 
that  an  electron  concentration  of  109  e/cm3  would  yield  about  10  db  (one¬ 
way)  or  20  db  (two-way)  absorption  over  a  20  meter  path.  Equation  I -4a 
would  require  a  40  meter  path  or  twice  the  electron  concentraton .  Simi¬ 
larly,  Equation  1-5  indicates  that  10l°  e/cm3  will  produce  high  levels  of 
absorption  at  a  frequency  of  95  Ghz.  One  thus  concludes  that  if  the  elec¬ 
tron  concentration  exceeds  109,  or  perhaps  1010  depending  on  geometry  and 
operating  frequency,  absorption  will  be  the  dominant  effect.  Thus  refrac¬ 
tion  is  only  of  interest  when  the  electron  concentration  is  below  this 
range. 

ELECTRON  CONCENTRATION 

To  maintain  an  electron  concentration  of  109  e/cm3  by  thermal 
collisions  alone  requires  a  local  temperature,  T,  of  about  2500*  K  for 
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sea  level  conditions.  However  the  ionization  generated  by  the  neutrons 
and  delayed  gamma  rays  from  the  fission  debris  can  maintain  this  concen¬ 
tration  to  a  much  lower  temperature.  The  gamma  source  is  in  general  the 
more  significant,  although  at  early  times  and  close-in  the  neutrons  can  be 
of  equal  signif icance.  Reference  1  gives  the  gamma  ray  ionization  source 
as 


where 


2  x  1019  Wp  p  -u/pdr 
y  4irR2(l  +  t)1 


ion  pairs  cm-3  sec-1  (1-6) 


Wp  =  fission  yield  (MT) 

p  =  air  density  at  the  field  point  (gm  cm"3) 

R  =  range  from  the  source  to  the  field  point  (km) 
t  =  time  after  detonation  (sec)  . 
p  =  mass  absorption  coefficient  (cm2  gm-1) 


In  the  following  we  are  interested  in  the  shocked  region  as  the  shock 
becomes  transparent.  This  occurs  at  a  shock  temperature  which  is  only 
weakly  dependent  on  yield,  thus  the  shock  radius  to  be  used  in  calculating 
q  scales  approximately  as  the  cube  root  of  the  yield.  The  fraction  of  the 
total  yield  which  is  fission  is  of  course  an  unknown,  but  1/2  .is  a  reason¬ 
able  nominal  value  for  megaton  class  yields.  This  fraction  tends  to 
be  larger  for  small  yields.  The  intervening  absorption  (e”w^p£*r) 
increases  with  yield  causing  a  decrease  in  q  as  the  yield  increases.  As  a 
result  of  these  various  factors,  q/p  in  the  region  of  interest  will  vary 
less  rapidly  than  the  cube  root  of  total  yield.  Furthermore  the  predic¬ 
tion  of  and  absorption  will  be  shown  to  vary  as  the  square  root  of  q. 
Thus  the  final  conclusion  is  only  weakly  dependent  on  the  inputs  chosen 
for  Equation  6.  To  represent  a  nominal  one-megaton  surface  burst  we  will 

choose  VL  =  i  MT,  R  =  1  km,  t  =  .3  sec  and  e-lJ^pdr  »  i  ,  and  obtain 
t  2  2 


q^  -  3  x  10 17  p  . 


To  approximately  account  for  neutrons  we  will  double  this  and  use 


q  «  6  x  1017  p , 


ion  pairs  cm-3  sec-1  .  (1-7) 


The  quasi-equilibrium  solution  of  the  rate  equations  for  the 
electron  concentration,  as  given  in  Reference  1,  is* 


where 


N  .  /!  'jlL  d 

6  /TJ5  ♦  0  *  A 


Aa^  +  Do^ 

a  =  -  , 

A  +  D 


cm 


-3 


cm3  sec-1 


(1-8) 


(1-9) 


-  ion-ion  recombination  coefficient 

-  3  x  10-8  +  6  x  10“6  _L_ .  ,  cm3  sec-1 

Tz 


(MO) 


ad  =  electron-ion  recombination  coefficient 


9  x  10-5 


cm3  sec-1 


(1-11) 


A  =  electron  attachment  rate 


-  9.7  x  103  exp{-  £22)  +  0.9  , 


sec 


-l 


(1-12) 


*  In  these  equations,  p  is  the  pressure  in  dynes  per  cm2  and  T  is  the 
absolute  temeprature  in  degrees  Kelvin. 
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D  =  electron  detachment  rate  which  is  the  sum  of  collisional 
detachment,  Dc,  resulting  from  high  temperatures  and  photo 
detachment  Dp  caused  by  energetic  photons  emitted  by  the 
fireball.  The  collisional  detachment  coefficient  is  given 
in  Reference  1  as 


D 

c 


2.4  x  10 


4  P 


—  exp(- 

/T 


5590 


)  +  2.1  p  /F  exp(- 


4990' 


sec 


_  l 


(1-13) 


The  photo  detachment  coefficient  is  given  in  Reference  2  as 

D  =  (^t)2  1.36  x  10-16  Tc5-4  ,  sec*1  (1-14) 

where  Tp  is  the  effective  radiating  temperature  ( °K)  and 
(^f.)  is  the  ratio  of  the  fireball  radius  to  the  range  to 
the  point  of  interest  (which  we  shall  take  to  be  unity.) 

To  compute  values  for  the  above  reaction  rates  and  estimate  the 
relative  importance  of  the  various  terms  we  will  use  ambient  conditions 
corresponding  to  surface  values  for  the  mid  United  States,  i.e.,  an  alti¬ 
tude  of  about  4000  ft,  where  the  ambient  density  is  about  1.1  x  10"3 
gw /cm3  and  the  ambient  pressure  is  9  x  105  dynes/cm2.  Earlier  studies 
have  shown  that  temperatures  greater  than  about  800°  K  lead  to  high  levels 
of  absorption.  We  will  use  a  shock  temperature  of  1000’  K  and  the 
corresponding  values  of  shock  overpressure  and  density.  These  are  an 
overpressure  ratio,  Ap/p,  of  15  and  a  density  compression  ratio  of  4.4. 
Then  the  pressure  at  the  shock  front  is  1.5  x  107  dyne/cm2.  This  yields 


A  =  9.7  x  103 


^•5  V?7!2  exP(- S2L 

(WP  1000 


)  +  0.9 


(1.5  xlO7)2 

1FF“ 


=  1.2  x  109  +  2  x  108  =  1.4  x  109  ,  sec-1 
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(Note  that  the  first  term  dominates  for  these  conditions) 


a  -  =  3  x  10‘8  +  6  x  10'6  . JL-12- 

1  (103)2-5 

=  3  x  10-8  +  2.85  x  10-6  =  2.9  x  10'6  ,  cm3/sec 

(Note  that  the  second  term  dominates) 

a  ,  -  —  —  1° —  =  9  x  10-8  ,  cm3/sec 
d  103 

D  -  2.4  x  104  1 exp(-  gg°)  ♦ 

c  /IF  1000 

2.1(1. 5  x  107)  /IF  exp( -  1£^£) 

1000 

=  4.25  x  107  +  6.8  x  106  =  4.9  x  107  ,  sec'1 

D  »  1.36  x  10-18  (10,500)5-4  -  106  ,  sec-1 
P 

The  last  two  equations  show  photo  detachment  can  be  ignored  (within  the 
shocked  zone) . 

Substituting  the  above  values  into  Equation  1-8  yields  an 

electron  concentration  at  the  above  specified  shock  condition  of  N  =  1.1 

0 

x  109  e/cm3,  showing  we  are  in  the  range  of  conditions  that  are  of 
interest. 

TEMPERATURE  DEPENDENCE  OF  ABSORPTION 

We  may  obtain  the  dependence  of  Ng  and  thus  absorption  on  local 
temperature  and  pressure  by  noting  the  dominant  terms  in  these  equations. 
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In  equation  1-9  for  a,  Ac^  dominates  the  numerator  and  A  dominates  the 
denominator,  i.e., 

Aaf  =  (1.4  x  109)  (2.9  x  10-6)  M  x  103 
»  Dad  -  (5  x  107)  (9  x  10-8)  =  5 

and  A  =  1.4  x  109  »  D  =  5  x  107 . 

Thus  Equation  1-9  can  be  closely  approximated  by 

a  -  a.  •  3  x  10-6  ( I -10a) 

We  may  also  simplify  Equation  (1-8)  by  noting  that 

=  [6xl017  (4.4)  l.lxlO-3  (3xl0-6)]l/2  =  9x10-  «  0 

Thus  Equation  (1-8)  becomes 

N«£^,cm-3  (1-15) 

e  A 

Now  by  using  only  the  dominant  terms  in  D,  A,  and  oi 
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we  may  us  the  gas  law  to  simplify  the  square-root  term  i.e.. 


-B  =  T  =  2.9  x  106  T 
p  M 


and  that  term  becomes 

[  ]i/*  .  2  *  108  T3/4  . 


The  resulting  approximate  formula  yields  predictions  that  are  a  few 
percent  higher  than  those  obtain  using  all  terms.  To  obtain  agreement  at 
1000*  K  we  will  use 

t3 .25  aQQn 

N  =  4.3  x  108  -  exp(~  - )  ,  cm'3  (1-16) 

e  P  T 


=  4.3  x  108  fill 
P 


where 


F(T)  =  T3-25  exp(-  49.9.°)  . 


(1-17) 


For  convenience  this  function  is  plotted  in  Figure  1-1. 


Equation  1-16  may  be  substituted  into  the  absorption  equations 
to  relate  absorption  directly  with  temperature,  pressure  and  frequency. 
From  Equation  (1-1),  (1-2)  and  (1-16) 


Ab  = 


8.5  x  1013  g  F(T) 

7.3  x  108  (g  p)2  +  (h  f)‘ 


db/m 


(1-18) 


when  u)  <  v,  then  this  reduces  to 
Ab  =  2  x  105  lill  , 


111  «  V  » 
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db/m 


(1-19) 


and  in  the  other  limit  of  <u  >  v,  then 

Ab  =  8.5  x  1013  lijl  ,  w  «  v  ,  db/m  (1-20) 

MAXIMUM  TEMPERATURE  FOR  REFRACTION 

In  order  to  specify  the  temperatures  above  which  absorption 
dominates  we  must  make  worst  case  assumptions  about  the  path  and  the 
system.  The  path  length  will  be  shortest  for  a  small  yield  -  but  even  the 
shock  of  a  5  kt  burst  is  about  250  meters  in  radius  when  it  becomes 
transparent.  Perhaps  one-tenth  of  this  is  a  minimum  path  length.  If  we 
also  choose  20  db  two-way  loss  as  a  system  limiting  factor,  then  for 
typical  radar  frequencies  we  will  use  Ab  *  -I  db/m.  By  inverting  Equation 
19  and  using  a  shock  pressure  of  about  15  atmospheres  we  obtain  F(T)  -  5  x 
108 .  Figure  1-1  shows  this  corresponds  to  a  temperature  of  1400°  K.  Note 
that  this  high  shock  pressure  also  implies  an  increase  in  the  applicable 
frequency  range  of  Equation  1-19,  since  v  is  proportional  to  p. 

The  above  temperature  is  higher  than  previously  suggested  as  an 
upper  limit  primarily  because  the  path  length  chosen  is  quite  small.  When 
the  path  length  through  the  shock  region  is  chosen  as  100  meters,  then 
this  temperature  drops  to  1130’  K.  We  have  also  assumed  the  absorption  is 
uniform  along  the  path  in  the  shock.  Since  it  is  not  -  i.e.,  it  depends 
on  p~2,  then  the  above  temperature  is  an  overestimate.  For  convenience  in 
setting  up  the  refraction  calculation  we  will  somewhat  arbitrarily  use 
1200*  K  as  our  upper  cut  off. 

There  are  several  assumptions  in  the  above  discussion  that 
should  not  strongly  limit  the  more  general  applicability  of  the  result. 

For  example,  the  ranges  considered  and  the  fission  yield  and  the  doubling 
to  account  for  the  neutrons  are  all  consistent  with  a  nominal  1  MT  near¬ 
surface  burst.  Shifting  to  a  nominal  small  yield  would  cause  several 
nearly  compensating  changes  in  the  numbers,  but  the  resulting  shock 
temperature  to  give  specified  db  levels  would  change  only  slightly. 
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APPENDIX  II 
SUBROUTINE  LISTIN6S 


pWUGWAM  HLAST 

0plVfcK  FP*  TEST  Of  ntfWACTTOM  EKRpM 
INC.  I  Ul-'E 

I  (HiJCA|  CONTINUE 

pP{  N  (  UN  1  T  «?,  N  A^F ;»  •  USF.*  1  jWfcFWACT  „Ot(T  •  .STATU  S«'Nf»< 
CONTINUE  ■  .FALSE. 

I  CASE  «  iA 

l'*  IC  ASF  «  ]  L  ASF  ♦  1 

1  Tf  ( p.  |j;«»»)ICaSE 

FC'WmAT(‘  INPUT  F  ON  C  ASF  NUMHfH  '#!•)) 

TALL  l»t’T  INPUT  (CONTINUE.) 

FALL  PFf  PACT  (UImF  TA) 

W«l  Tf.  (  ?,  1  hNujpThE  t  a 

lv1»>p  fUnKAIC*  Uf f  R  ACTIVE  f-KMOfi  IS  '  »1pFl-J.b/) 

AMllHtk  SIT  OF  INPUT  OAT  A 
IF  (C'OjTTMIE)  NO  TO  IP 

FNf'  OF  IMPi'T  * '  A  T  A 
CLOSE  Cun! T  «? ) 
f  NO 


oonc*  o  o  o  o  o  o  o  n  o  n  nn  n 


SllHROUTlMt  (,FT  lNPDT(MOrtt  DAT*) 


BEADS  INPUT  F  WDM  WEFWACT.DaT 

TRANSEFBS  THE  I)  A  )  A  VJA  THt  COMMON  MLOO  f;IVFN.CMN 

LOGICAL  F 00  A|.  ,  MDkf  l>A  T  A 
TNCLUDf  *1,1  VFN.C^’V* 

D  I  MENS  TON  IHOLCHM)  ,  I*(2),  i9H()A(‘>),lTfcMPACM,T'"ATE«(<>), 
1  TWHSC4), I«TSf 4), JThETAI ( 7 ) , I CHMMfc My ( 3) , 

?  TLMDC4) , IT lMf (5) 

DATA  J>/1,1hw/, 

1  lPHOA/4, \H(i,  1HM,  1MO,  IhA/, 

2  ITEMPA/5, IhT, IMS, 1hm# ^HP, IHA/, 

3  UATLP/h,  ! MM,  IHA,  1hT,iHF,1m«,, 

4  lPbS/3, 1 HW, 1 Hh, 1 hs/, 

*»  IHTS/3,  ImW,  luT,  IMS/, 

6  iTMFTAJ/f),  l*iT,  !HH,  ihL#  ImT,  1Ha,  |H1/, 

7  ICOMMFNT/?,  IK,  1M?/# 

M  1FMO/J, IMF, !HN, 1MH/, 

9  iTlMfc /4, tHT, l  HI, |HM# ImE/ 

data  whoa /i  , ??bF« J/,  tf  mpa/?wp,  /, *  a tew/o ,  m  / 


IF  (MUWtfATA)  GO  TO  >  ” 

pPENtlJNlT«t  ,  NAMf«  *MSfWl  *Wf  FH  ACT.  0*T  •  ,  ST  A  TUS»  '  01.  D  *  ,WEAI)QNLV) 

C,ET  DOF  |  INF  OF  INPUT 

10  BEAUT  1  ,  lHbD,F_ND«b  H  )  IhuL 
1  P50  FOWMAT I8da 1 ) 


SKIP  THfc  L  f*>r  IF  IT  IS  A  Cf)MMfcNT 

CALI  MOLF'Ji  Cl  hoc  ,  I  comment,  f  DUAL  ,  ,F  ALSt  .  ,  NS*  IP  ) 
JF  (EUUAt  )  GO  To  M 
C 

c 

CALL  MOL FUL ( I  MOL# I  WHO A, f  uu#L,.T«Uf.»NSK IP) 

Jf  (EUUAl  )  Ml  To  1  M 

C 

CALL  HOLFUI  CIMU  ,IHMF>4,f  JoA|  ,  ,  T«uF  ,  ,  NSK I P  ) 

TP  CfcuHAl. )  <;n  To  l?o 
C 

CALL  MOLFUl  (  I  Mill  ,  I-ATi  w,F  Dial,  ,TWoF,,NSKIP) 
fF  ( F u 1 1 A |  )  t;n  To  f  3  i 


44 


nnn  onft  oon  o  r>  o 


C 

CALL  MOLgUl.  (IHOL  ,  IRrtS.EuUAi  ,  .  TRUE ^  ,  NSKJ  P ) 

IF  (EQUAL)  60  TP  17* 

C 

r ALL  MOLFUL ( IhMl, I»TS*FuM*t , . T RUE . , NS* I P ) 

IF  (EQUAL)  GO  To  1 
C 

TALL  H0LFQL<1HU|  , I  The T k\ , t QUaL , , T  rUE , , NSK l P ) 
IF  ( k  uu A|  )  (,0  TO 
C 

TALL  MOLEUL  (  Th'IL  ,  I  T  lMh#tUUAt  #  ,  TRUp  .  ,  *Sk  IP) 

IF  (E UnAl  )  GO  TO  ?♦)■* 

c 

CALL  MOLFUl  (IhUL,  !•>, EQUAL,  . Thue.#m3KIP) 

IF  (Eui)AL)  GO  To  1 

end  of  data  set? 

fall  HOLFUI.  ( 1HUL  ,  IENL),EOUA|  ,  .FALSf.,NSkIP) 

If  (EQUAL)  GO  TO  5*0 


FRROR  -  UNPEtl'GQlZAHLF.  INPpT 

RRITF  l  2#  1  *  0/  )  1  mol 

|06t»  FORMAT  (  »  UNRFCUGNWAHl  E  INPUT  '  (S»,H.'A|| 

STOP 

READ  w  (VIF.LO) 

100  BACKSPACE  1 

RE  AD (  l ,  1 P70 )  * 

1070  FORMAT (<NSKIP>X,£PV,V) 

TF  ("«l-  T  .  |  ,*f  «*.1)  wRI  TE  (2,  1  n75) 

1075  FORMAff*  R  ARN  I NG I  INPUT  YIELD  IS  |E3S  THAN  I  TON') 

I F  (w,lF.0..O  Then 

mRITF(2,TU 7b)  * 

1076  FORMA  T  ( *  ERROR!  VltLO  ■  *,PFli«.4,»  IS  NFGaT/Vf  OR  ZERO.  PROGRAM 
1  EXITS.') 

STOP 
FNOIF 
GO  TO  lU 

READ  HHOA  (DENSITY) 

||0  BACKSPACE  1 

RE  AU( I , I07P)RHUA 

TF  (RMfTA.LT.  I.  *f-6)  tPfti)) 

1 080  FORMAT  ( *  «APNING|  INPUT  OInSITy  1*  Lf  SS  THAR  IF. 6, 

1HLAST  model  IS  SUSPECT,  •) 

60  TO  10 
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c 

c 

c 


PEAO  1 F  MPA  (AMhJtNT  TEMP.) 


120  hAOSPaCF  » 

MEAD(  )  ,  )i<7P)  TFmPA 
T F  (ItMPA.Ll  .lVHi)  THEN 
w»lTF.t?,1l  #**S  > 

1HHS  EOMMAT<*  M  A  P  M 1 0  (» I  INPUT  AMrIFNT  TpMP  IS  VEWY  SMALL  0« 

jnegative.  nom i n al  value  he  2*bk  used.*) 

TEmPA«?8«. 

FNDIF 

GO  TO  to 

c 

C  Pt AO  nATEP 

c 

tJ*  BACKSPACE  t 

PtAnn  ,  |tl7«)  **■  A  T  k  M 

TF  (mATFm,LT,w,0)  THCn 
wMITEO,  lua.i) 

1H9-A  EOMMATf*  PAWNING:  INPUT  VLAUE  OF  wATEM  VaPOK  IS  NEGATIVE, 
JOHV  A ] H  TS  USE  P,  1  ) 
wA  Tk  P»P  ,  w 

END  IE 
GO  TO  10 

c 

C  BEAD  KHS  (MANGE  EMOM  SOUPLF  TO  HUpST) 

c 

17*  PACKSPACF  1 

ME  AO ( 1 ,  I  o  7^ )  MBS 
TE  (MHS.I  T  , O . v»  T  THEN 
WMITE12,  l.’O?)  MhS 

1*92  EOPMATC*  EPMOWi  MANGE  EPUM  SuUMCE  TO  Bl'MST  ■  •  ,  PF.  J  tf,  4 ,  *  JS 
1  NEGATIVE.  PWOGMAM  t*ITS,') 

STOP 
ENOIE 
GO  TO  1* 

C 

C  PEAD  MTS  (PANGh  FMiim  SOUPCE  TO  T*pG£T) 

C 

18*  BACKSPACE  1 

PEADM.  |*7P>  MTS 
I E  (HTs.IT.k.a)  Then 
mKITE <?, 109b)  HT3 

JM95  E  OMM  A  T  (  *  EM  MOM  I  MANGE  F  MOM  SO'JHCE  TO  TaMGET  ■*,PElW,4,«  IS 
1  NEGATIVE,  PMOGPAh  EXITS.') 

STOP 
END  IE 
GO  TO  |0 


c 

C  read  T  Hf.  T  A  t  (ANGLE  HFTwFFN  VFCTORs  WMS  AND  RTS) 

C 

19*  haCKSPaCE  1 

REAU(1,1*7m)  T  *C  T  A  1 
IF  CTHFTA.LT. *.*)  THfN 
HRITECP,?^  )  ThL  T  A l 

2HMO  FO«mAT(«  ERROR!  THfctAl  a  '  ,  PE  1  * ,  4  ,  •  IS  NEGAT  IVE,  PROGRAM  EXITS,’) 
STOP 
FNOIF 

GO  TO  {to 

C 

C  READ  T  trip 

C 

2t4*  BAOSPACF  1 

HE  AO(  1 ,  Ut7fl)  T  I  ME 
TF  CTImf.LF .*.*)  then 
*RITE<?,?wh5)  T  i  i^F 

2*t'R  F  OR M  A  T  (  •  ERRORl  TIhl  *,Pfctei.4,'  tR  NEGATIVE.  PROGRAM  EXITS') 

STOP 
END  IF 

r,n  to  i vj 

c 

c  FRO  STATF^FfJT 

C 

50L  fOrtl  IMtfc. 

MC1REOATA  a  .TRUE, 
no  TO  AS* 

c 

c  FNI)  OF  FTLF 

c 

6WP  rOUTIMlh 

MORkUATA  ■  . F  A I  S f . 
rL  USE  (OUT  T  a  |  ) 

GO  TO 

C 

c 

t>J>*  rR  1  T  F  ( ?0  oiuv*) 

9*&m  FORMAT  f  *  * ) 

wRT  TE(?,0  ffh) 

9*fl!>  F  OHM  A  T  f  '  Input  data  •> 


KtUTt  l?,  0171  )  * 

9*71  FOHMAlf*  w  CVFItD)  . . . SlPflJ.h) 

**«I  TF  (?, OH72)  RHUA 

9P72  FORMAIf*  R»u-a  (dEnsTTT)  . MPEI3.A) 

h»ITT  (?,Q*'7J)  Tl  mpa 

9*73  FORMAtf  TFmpa  (AmhIIOT  TEmP.)  . '.1PEI3.6) 

►  RITE ( 9, Ot  74 )  MATER 
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9*74  FORMATf  water  . . . 1  Pt  1  3 , 0  ) 

wRITE(2,9U78)  RH$ 

9H7ft  FORMAT  (  •  RMS  (hijRsT. SOURCE  t> I  S T  ANf  E  )  1  Rfe  1 3 . t> 7 

wRITE  (2*0079)  RTS 

9«79  FURMATf*  WTS  (TARf.fcT-SOliRtF  l)IST*MCfc)  .  .  .  •  »  IPE  t  3.6) 
wP1TE(?,9*sm*)  ThETAI 

9M8M  FORMAT  (  «  THETAI  . . . 1  #  IPE- 13,*>) 

wRITET2,9j*«n  Time 

9M81  FORMAT^  time  . . . * »  I  PE  1 3 . o ) 

WRITE 

C 

PE  TURN 
END 

C 

c 

SUBROOT  I Nh  HOI.FJJI.  (I  ONE  ,  J  T  wn  ,  ANS  ,  l  p NG  T H  ,  I SK  I P  ) 

MODIFICATION  UF  mRCSIm  SUrROUT  IMp  ;  TESTS  F()K  EUUIVALEHCfc 
HETwtEN  ARRAYS  I  ONE  AND  IT^OT  IF  LENGTH  IS  .TRUE,,  IT 
RETURNS  The  NUMBER  of  CHARACTERS  FROM  THE  Hp  G I NN I NG  UF  THE 
RFCORO  UNTIL  AU  *  CHARACTER  IS  WfaChEO, 

I UCICAl  LENGTH,  ans 
DIMENSION  IONE  f  I  1,  ITwlifl  ) 

OATA  1HLANK/JH  / . I E  UU/ 1 H • / 

ANS«.F  ALSt. 

TSKlP-v' 

J«P 

NtlTmjf  11*1 
no  i*  io.n 
2  J  ■  J  ♦  1 

IF(IONEU).EU.IHLAHK)  GU  Tn  ? 

TF  (  IONF  (  J)  ,NE,  I  TRU<  I  )  )  r.U  TU  UM 
1«  CONTINUE 
ANS« , True , 

TF(LfRRTH)  GO  TO  IS 

RETURN 

15  CONTINUE 
K«J 

DO  2*  L«*»8v 
.!■  JA  1 

TFC  I  0  N  p  (.])  .  T  N.IfUll)  GU  TO  AO 
2B  CONTINUE 
3*  TS*iP»J 

IP* A  RETURN 

end 


w  ^ 


SUBROUTINE  REFRACT CUTHETA) 

C 

INCLUPF  MiTVEN.CMN» 

INCLUDE  *wFRT.CMN* 

DATA  RM07/1  .3/,TtMPZ/;>8R,/ 

HALFPI«ASIN(  )  .vi) 

PI«2*HALFPT 

C 

nTHETA*,?..) 

wJ*ft»*RHOZ/RHnA*TFR,PZ/TFMPA)#*nt/i.) 

STIME«TIme/w3 
CALI  DFNS1T*(STI*E) 

C 

C  1 1  ST  Flip  t On l /AT  I  UN 

IF  (TEMPK.LE.  12h*. )  GU  TO  tB 
wR  I  TE  (  2,  fMM-i  ) 

8H0V)  FORMAT  r»  SHOCK  STRUNG  ENOUGH  To  pAUSt  ABSORPTION,  REFRACTION 
1  IGNORED •  ) 

GO  TO  8*  U 

c 

C  «1  is  tmf  Range  from  source  Tn  point  of  closest  approach 

|H  P1*PHS*C0SCTHF  Tat  ) 

SRA0«PPAD*W J 

RTB«S'J9Tf  RHS**p4RTS**?«2*RPS*RTS*rl,Sf  THETAJ  )  T 
RHIN.hHS*SlN(THF  T  A  1  ) 

L  AST  !Mp-v» 

c 

C  IF  SHOCK  HAS  NUT  RE  A  CDF  D  LINE  OF  SIGHT,  r>n  NOTHING. 

IF  (RHIN.GF  ,SRAO)  GO  To  Air. 

C  IF  SHOCK  HAS  NOT  CROSSED  LINE  OF  SIGHT,  DO  NOTHING. 

TF  (  (RTH.GF  ,8RAf»)  .AND.  TRTS'.L  T.RJ  )  >  GO  TO  4" 

TF  ((RRS.GF  .SWADT.ANO.fRt,(  T.*))  r,0  To  4in 
GO  TO  SO 
4«  CONTINUE 

WRITE  (2,801M 

8WJW  FORMAT (•  SHOCK  HAS  NOT  REApMfcU  LInF  OF  SIGHT.*) 

GO  To  «ua 

C 

5P  CONTINUE 

SIPSID«(PHTn/SRAO)  I 

PSIP*ASIn(S!PSTf ) 

FTAM«1 .♦«.P2*PMOA* ( 1 .*4.BF ♦3*WATEp/TF  MPA) 

C  DEFINE  etas, CUNMAX, AND  CONST 

STRT»AMINMSPAO,AHAX1{WhS,PT«)> 

STRTS»STPT/w3 

CALL  ETA (STRTS,RHUT , TEMPI ,pT AS) 

JME  »AM  f  N 1  C SPAN, AMIN) (RRS,Wxh)) 

8ML  S«SMl,  /M 
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CALL  t  T  4(SMLS,RH')I  ,  TfMHi  , E  T  AM ) 
rONM*x*SML  $*E  T  4m 
TF  (SlRT .tO.SRAo)  THE.N 

CONST»PraO*E.TA,  *STPSI  4 
FLSF 

IF  f Sml.FO.WmS)  ThMi 

rowsT»roMM4x*siNf the ta  j ) 
tLSF 

SFPS«WHS/w3 

fall  ETa(swks,Rhoi ,  Tf  MR  J  #E1  iSRHS) 

rOMbt«RHHs*t  T AS^MS#StH( THf  r4  t  ) 

ENDTF 
FND  I  E 

6«  LOOPC*P 
C 
r. 

IF  (SPAD.RT.AMlNt <RbS,RTR))  MJ  TU  9 <4 

f  U  m  P  |  F  T  f  TRAVERSAL  .  HoTm  TaRrFT  a  Kit)  SOURCE  ARE  OUTSIDE  SHOCK 
INTER  KATE  FROM  S^AO  TO  RM  I N 

PfcCALC  »hjm  MS  f  Mi;  f  T  A  aT  RHJN 
SRMIKsRMTN/RS 

CALL  ETA(SRmTN,rh.im,  T tMPM , F T ARM  I N > 

rmin«rmin*ftarmim/etavi 
7  5  DU*I«PPAI.) 

nUH?*KMl N/W  J 

FTANEw.ETAS 

CALI  InTFC.PT  foUMI , OUN2, CONrT, ETANpR, Ph I  ,.J  ) 

FTAPMIM.FT ANEW 
F  TAHAK.CONST /PR AD /COS (PH  I) 

nTHf  TAiCSRAO^SORTCETAHAHI-pMINl/f pHS ♦CDS ( T HE T A  1 ) ) 

tF  clasTthf  .eo.  l )  r.n  To  RHP 
SIPSI^«RPS/3RAD*SIN(THETA1 aOTMETA) 

PSTPBASINCSIPSI;'* » 

L  OOPC •LODPC ♦ \ 

CONSTk^CONST 
r0NST«ETAi'*PPAO*SIPSI)' 
rhjnbCUOST  m<.3/T  T  ARhIn 
IF  (CONST. RE. CON MAX)  Then 
CONSTsCONMa* 

RMlN*ri'NST***4/F  TARN  IN 
LASTiMt.l 
ENUTF 
60  TO  7b 
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c 

L  INTER  M*TF  *Lt  CASES 

C  TAMGfT  iiW  Si'UHCF  UP  HIltM  APF  iNSlDF  SHI'(-K 

yw  ThF  TA3aASlN(kTS/PTF>#SlN{  THp  T  4  1  )  1 
C  ASfN  WFTt/PWi  t/AU'FS  HtT«FF>  -hALEPI  A*W  MALE  PI 

f.  ChfCK  IF  TMfT43  SHOULD  #F  UPEaTEP  Than  HaLEPI 

T F  ((P|,|  E  ."T.'iP.CRTS.Lf .«m  GO  T n  9b 
TME  TA^«ASI^(khS/HIh*SIN{ THpTAl ) ) 

IF  (( ThF TAl ♦Thf T*?)  .1  T.HALFPt )  THrTA3aPI,THtTA3 
C 

9b  rl'NTlNilf 

10**  OliMl«AMjNnWKS,SwAI.J/H3 
|)UM?«  Ah  I  m  J  f  PT*3,  SPAO)  Z"3 
F  T  ANtwaE  T AS 
nPHin*«i  # 

FALL  InTFOPT  n»>D’|  ,l/NM?,rDNsT,tTAN(rW,  PHl,OpHlr>K) 

T F  (IPl.Lfc.  >I.DW.(PTS.LE,f<l))  GO  TO  ISO 

LINE  CF  SIGHT  CWl'SSFS  POINT  Dp  CLOSEST  APPMOACh 
it'll  I  T  i  ON  At.  I  o  T  F  GkATIliO 
OHM  1  mAMlMl  CPPS,»'  TM)  /*.% 
f)OH?»HM  I  n/w  3 
nPHlOKH* 1  . 

CALL  1  N  TFRb  I  I  |,‘,M  1  #  t'OH^,  (  ()H«jT  ,  t  T  AOpw  ,  Ph  I  vl  ,  ('PH  r  DA  vl ) 
PHJ  a  P  H  J  *^*Ph  I  it 
OPHIDA  (t/PM  f  OK  aPaHPH  JI.»KO 
pTAPHji\,»fT  AV'F  M 

C 

C 

lb  ’  ft  (PHS.tT.3PAD)  GO  To  j  7 b 
C  SDlikfF  IS  DuT  Slot  and  tAWGFT  is  INSIDE  ShOC* 

l  ■ SoK  T (KH3**?*3PAD**2»?*«Hs*SwAO*rnS( THF TA3.PHI  )  ) 
oThE  TA«  A  SI  N(SPA0/!.«SIN(  THE  TA  3-PoI  )  1-Theta  l 
IF  (LAST  THF  .FM,  |  )  GO  TO  lu'fi 
PS  I  E  a  T  HE  T  A  T-PH|  ♦  ThF  T  a  1  ♦(> T Up  T  A 
C 

C  l Tf W4T TON  l NOP 

i  oopc »i  ijoPf,  *  t 
CONST  ✓) a l n n S T 

rONSTaPNAO*E TAA*SIN(PST  1) 

PH  I  NaCONST  **3/l  TAPhJN 
TF  (CONST. GF‘. Cf*NMAX  )  T  HE  N 
CCNSTaCONM A* 
hHiN»C0NST**»3/fc  TAH^O' 

l  AST  I  Hi,  at 

t  NO  TE 
GO  TO  | 30 
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c 

c 

175  IF  (HTh.l  T  .SPAD)  &0  TI)  ?v*h 
C  source  IS  iMSli-fc  and  T*wgit  1*5  OUTSIDE  Shuck 

C»SUH  r(KTH«*?*S«A|)**2-?*KTH*SPAn#rtiS(  THFT  A  J.PHI  )  ) 
PSIM«ASIN(«TH/C*SIN( The  T A 3-PhI H-hALFF  I 
IF  CLASTIHE  .Eu. 1 )  GO  Tp  l«p 
C 

C  ITERATION  L OUb 

L  OOPC  »L OPPf ♦ 1 
CONST  4»COnST 

C ONS T  «PR au*E  T  Av.-bS  1  n  (  PS  J .; ) 

MMIN«CnNST*w.T/t  TAWM1U 
TF  (CONST. Gf .COMHAX  )  Tnr  N 
CONST  *CuNMAx 

WPlN»CONST*'t3/FTAKM|^ 

LAST  I Mk  *1 
E  NO  t  F 

go  ru  iww 

c 

lflfl  SRb8«PHS/w3 

CALL  ETA  (SHHS,  #T£  IP]  ,k  TARHS) 

OTHET  A*ASlN(  TSRAU*kT  A,t*slNfPSI^)  )  /(RHS*LTawhs)  )-TmETA1 
GO  10  tri-A 


SOllWCf  ANP  TARGET  tut  T  NS  lot  SHOCK 
2bb  IF  CL  ASTTMf  .FU.l  )  GO  To 
C 

C  lTfHiTTOM  LOUP 

c  OOpC  *1.  U0pc  *  1 

CONST  »t«CONST 

CONSTaf  ThETaj-PhI  J/hPMlPK  ♦cOnST  i 
PMIN«COHST*w.i/F  T  Ahr  I  n 
IF  (CONST. GF.  .comma*  ) 

CONSTbCOnha* 

Pm1n«c0NST*«3/LTARMIn 
l  A  ST  I Mt • 1 
ENUIF 
GO  TO  1'id 

c 

25b  sPPSbRhS/h3 

CALL  tTA(SRHS,RHOl ,  Tt»'Plf k  tAPMS) 

nTHE  TAiASIA'C const  /(«PS*F:TAwKS)  J-ThETai 

c 

8bM  RETURN 

FNP 
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SUhRlK'TlMfc  IMT£GwT  ( WM,Cfc  ,E  TANfcw#PMl  ,l>Pl>*  ) 


INPUT 

MS  ■  LIMIT  0*  lNlfcf.MATlON 
MM  a  LIMIT  Of  INTEGRATION 
CK  a  CONSTANT 

ETANfW  ■  INDEX  nF  PEFPACTlON  AT  THE  UPPE*  L  I M I T  UF  INTEGRATION 
OUTPUT 

PHI  ■  ANGLE  METwLEN  THt  VtCTuRS  MS  ANU  PM 

OPD*  ■  THE  DERIVATIVE  Of  PhI  wITH  PESPF.CT  TO  THE  CONSTANT 


INCLUDE  *«ERT.CmN» 

INfLUUF  'GIVEN. CMN* 

data  num/sm/, testl/k'.ui /# tests tt, a"2/ 

c 

f  T  Atf •£  T  AnE  w 
RlaAMAxJ  fKS,RM) 

M2«AM1N1 (PS#PM) 

I*p 

SUMl*^ 

SUM2«t« 

C  IF  IETAH.NE.E TAC) 

STFPalPl-M?) /FI dAT(NUM) 

RADI aP| 

IH  RAD2«AMAXl  IP?, MAD1.STEP) 

IF  (PAD2.tU.R2l  STEP»WADI-«A!)2 
»ADCa(RAf)?*RAP1 )/?. 

1*1*1 

TF  Cl.r,T.Nl|M*2l  go  TO  2** 

CALL  f.TA(MADf  ,X1#*2,F.TAC) 

C 

TEST»AH5(tTAii.FTAC)/(tT  AV-i  .) 

IF  (TfcST.LP .TESTL)  NO  To  2* 

C 

9TEP*STEP/2. 

GO  TO  lv) 

c 

20  Ci*CK/ET AL/PA01 

C2«AMIN1  (1  .,CK/FTAC/AMAXHP2,MAI>2) ) 

ni*Acus(cn*AC's(C2) 

9UMl«8llMt  ♦n | 

IF  (AMaXI  <Ct,C?).Gt.l.)  GO  TO 
IF  (DPdk.LT, M, )  GO  TO  2*, 

D2*l./SQPT(l ./Cl **2-1 .1-1 ,/SuRT(r./C2**2*l.) 
SUM?«SUM2*D? 

C 

25  IF  (RAD2.LF.R21  G'»  TO  J( 

£  T  AH»t T  AC 
RADJaK AD2 
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SUHHfiUTlwt  E.T*(SW#whu^,TEmPR#6TAK) 


INPUT 

SW  •  SCALE"  RANbf  (C*) 

OUTPUT 

WHOP  •  f)FNSr  TV  (UVNC/CW2) 

U.MPP  ,  Tf-  upt«ATUWE  (DfGME.ES  k) 

ETAR  a  IHiFx  ilT  REFRACTION 

CALCULATE  S  RHOP  A  No  Tt^Pp  FROM  1  k!  STANDARD,  ANl) 
THfcfM  FTAP  AT  SCALED  PANGf 

INCLUDE  'GIVEN. CMN* 

|  Nft  UOf  *Uf»'SlT  r.CPN* 

DATA  Rho;/ 1 .225E-3/ 

CALL  OENS(SR) 

WHOPakhua*! 1 .♦OoW/RuuZ) 

TtMPR«A«A*t ( TT  mpa,TEMPK*{ { 1 ,*OUM) /( l.*PDPK J 
F TAPal 1 .♦4#*E3*WATE«/TEMPR) 

RE  TURN 
F  NO 


r>r»r»ononr»oooor» 


$Uf»KUl)TIKj£  »If?PT  (£fcr,»KW#r,H()Nf  ,PpFS,  T£MP) 

OOAN-NICKLE  FiJi'AMU*  OF  ST*TE  OF  aIH  (SEhI-PhVSICAl  FIT) 
AS  F*TWAt  T»  0  FhOM  MoaC  PEP  H?  vFw»lr>N  OF  [  AmH 

INPUT 

FEE  ■  FNFPGV  (EHGS/Gm) 

«PW  •  OENSTTY  (UH/CM3) 

OUTPUT 

GPONf  •  GAHhA  .  | 

PWEJ>  •  PPESSUPk  (OVNf/C'1?) 

TEMP  ■  Tf«PtPATUKE  (f)EGHFES  *) 


IGO*  1 
GO  TU  1 

c 

C  THIS  t  NT W Y  HE  TOWNS  UNLV  f»Ni|NE 

FNTRY  AlP(FFF,MPP.GM0"n 
I G  0  ■  P 
C 

1  F«EEE*  l  ,F-1v\ 

BH0» Wkp / 1 

IF  t  AbStFIT  ,L  T  .*>.U)  60  TO  .1 

c 

IF  (E  I.GT.P.td)  GO  TO  2 
C 

FO»rt. a 

FON«EXP(.E /o.h)) 

wS*0,vi 

GO  TO  4 

c 

2  F0*EXP(«F/4,4  r>) 

FQN»P,<* 

wS« I 
GO  TO  4 
C 

3  nM*.975*(RH0)**.ph 

FE1  t‘)*i»t4?*Al06{WH0) 

F I ■ C  EE t-F)/uFI 
mS»I ,/CtxPT-F I  )*1 .) 

F0«t*P(.F/4.4M**S 
FON»E*P  (.£/*.*>  3)  *<  J  ,»wS) 

C 

4  pfcTA«!-«.U 

TF  (t.LE.1.4)  GO  TO  !S 
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F2»CE»4k).v)/S.^ 

JF  (ABS(£2).LT.5.^)  t»0  T 0  7 

IF  (F2.GT,**.f1)  GO  TO  O 

F  N  •  n  ,  .1 
wS«W,H 
GO  TO  * 

FN»txP(-F/?b.M 

I^S*  1  .  H 
G"  TO  R 

nF?«4,n*pMo**n,M»b 

FE2»45,t»*«HO**u.*157 

FP«(E«F£?)/bF2 

wS«l.»/(EXP(.t2)*l.M) 

FN«EXP(-F/?S.5)*"S 

F3*(E«16w,m)/6.H 
HETA»BFT  A  *  (  l  *i4b*WS 

Ft**),«1 

IF  (E3.GT.f*5.*M  )  F  t »  1  .**/ (FXP  C -E3)  ♦  I  .  d ) 

r,HONF*{  ,161  ♦.?5h*F0*,28*F0N*.l  37*f^*.H5*FF  )*«ho**HETA 

IF  (IUU.EU.tf)  GO  TO  1*1 


PBES«GMOMt  *EfE*PWP 


TEHPtH*TO«E 


»HOLN»ALDG(PH(l) 

F«(b.b82549F.*i5*«HQLN-2.7l4  34F.id3>*RH0LN*9.7pt»0I 
J3«7fl«*E«!S*BHULN-t>.4lRb73E-4i*KMOlN*2.64bf  •? 
h«(  3,9231  23E-0  7*«HULN*5,97i  54 9F»t»i*HH0LN-9,?lt»Pb 
CONI »34fiW,P*GMONF*E 
eON2«(H*F*G)*F*F 
RE  TA»(E»3,M) /n,6b 

IF  (HLTA.GT,  Ji1.il)  GU  TO  V 

TEWP«CnNi/(C0N7*f  J./i*COfi7)/(FXP(HpT4)*l  ,p)  ) 

GO  TO  )k) 

T£MP*C0N1 /CON2 

RE  TOWN 
FNO 


«i  «  ■TtWifi-'Ar  f' 


o  n  o  o  r»  o  o  r»  ,-innnnonnnnnnnnnn 


SUhWOUTlNt  OfciSlTY  ( T  ) 


PKOV  1  Up S  MAST  PARAMETERS  AT  INPUt  TIMf  •T"  rAvSED  ON 
NIJtLtAR  MAST  STANDARD  -  1  kT,  Ap*L-TR-73-S5>  REV  APRIL  1975 

Initial  rminjTjuNS  assumed  to  he  stp,  i.e. 

AMHIFNT  PRESSURE  *  1.0t2sE6  DYNF/CM2 
AMHIFNT  DENSITY  a  1.325F-3  GM/rM3 

OUTPUT  VIA  COMMON  HLOCK  On  initiae  CA| L  AT  TIME  T 
PR  A  D  a  SHUCK  FRONT  RADIUS  f  CM  ) 

OPPk  a  Snjf>  FRONT  UVE«-PRESSURF  (OYNE/CMJ) 

0 0 P K  a  SHOCK  FRONT  OVER-nENSITY  (GM/fM3) 

TEMPKa  SHOCK  FRONT  Tf  Mpfr PATURE  fK)  -FROM  EQUATION  OF  STATE 
OUTPUT  ON  simSEOOFUT  CALLS  AT  SAmF  TIME  VIA  ENTRY  DENS 
OUR  a  OVER. UfNSITY  AT  RaOIUS  R  (GM/cm3) 

INCLUDE  »WFRT.CMN» 

data  T  A3 /-?*.?•>/,  TPM.  1/0.371  /.TP*p2/W,79/, 

1  AC/3,  IMF  I  «/,  AO/!  .MME14/, A$TAp/9.«F 9/, R STAR/4. 454E4/, 

?  RHOZ/1  .225E-W, 

3  H/O.H329I /, C/»!,RHh/#rZ/33S9T,M/,H7/849o,0/, 

4  HR/tl. 03499/, CC/-1 ,06H/, 

5  T  TUI .0/3.0/ 

TF  n.EW.TTOLU)  GO  TO  9 
f)E T ERM I NF  WFPR2  (RADIUS) 

1  IF  (T,lF.u.;i)  GO  *0  14 

calculate  Radius  Explicitly. 

early  time  form. 

OE  Tf  RM I Nf  RJP 

RZPa( 1 ,r-h*T**C> *(C/aT*hZ) 

C 

TF  (T.LT.P.2H)  wFPR2a242IU.*T**TPwRI*Ct,A(l,?3AT*P,l83)* 

1 C I.P-EAP( TA3«Ta#1 P*R2) ) ) 

TF  (T  ,1  T.w.n  GO  TO  2 

c  late  time  form, 

ALT  BALOG (T ) 

ALFTaALT*3.39|M? 

PN£HaR/P*| .U45E4A2, 15L3bAI Ta(9.E3/ALFT.0,8E3)/ALFT 

c  intermediate  time  interpolation. 

IF  (T.|  T.0.2H)  RUE^alRNEpaf T*H.| )*wFP«?*(W.2P-T))/M,18 
l»f  PP2*RNEh 
?  RewFPR? 
pRADaR 
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o  o  r*  r>  o  o  o  o  non  n  n n 


PETtNHIiME  OPPk 

3 

PT10«2.24f>17F-8*W 

CP«$U*T (  AlPGfPT  I0*3,^*FXPT.(S'JP1  IpT  X O )  /J.tf)  )  )  ) 

OPPK*  f  TAC^PAAMJAPKAAbTAP/fF  )*WP 

DETERMINE  POP* 

4  pp»oppk 
BTlU«t'P/l  ,.'1250 
P«0P*1.->12,SE»5 
6M0fc|«*.4 
r,  a  h  m  a  ■  i  #  a 

GAhkaigamha 
rtn  5 

PMC  I  »P  Mi)/ A  (  <P.*r,4**RA*  (GAM^A*  1  .  >*RTlO)/(2.*G4MRAAn,4MNA-l  ,  )*f?T  IU)  ) 
ee«p/(c;mcmf*pmoi  j 

CAUL  AJPPT  ((■  t  ,PHP1  ,pPF  s,  TtMPK) 

RAmRAC'«GA*PA 

GA"RA«2..*«r;Kt)Uf  /(  2.*>*Gm(, ht  ♦  I  .»')♦! ’.P 
TP  (  AHStGAMPA.r.AMPAO)  .1  T.  I  >»■ -M  R>-»  TP  * 

5  CUNTIMIE. 

6  POPK  »Km(j  1  *KM(V 

PE  TERr  JNf.  w/p 

8  PZPiM.ffl 

TP  (T.LT„i!.?6M  PZP«?,*bMt4*T**n.-,t)S 
TP  (T, RE. p.288)  KZl)»(l.P-MH*T**cr>*(CZ*T*8Z)*?PP,u 
Pf  TUku 

FMTWY  OPM5(MAi>) 

P«PAO 

PETFKE'TNP  OUP 

9  rPk  »PR  Af) 

IP  (T«L !.»•?)  UP  Tu  11 

C 

pPr  «WPk *  1 .05 
1  ,05 

TP  (  T «f  U.  T  T PI. !) )  GO  TO  1  *i 
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c 

P/aWZo 

pmnbUZ.9,71 bit 1*T**.I21|5 
wZ*PZ*t  .F-* 

PMNaW^N a l .f  .P 

0DMNa»t»,*1*tmPKA?.PF -5*T#*  (.1  .bvl?b) 

PNE&aNJ.pMN 
PPL  SsNPK.P/ 

HNPbNPk.HM^ 
it  M*U(»PK  /MPLS 

hln»uupk.aln*ppk 

ODMNLN«AL^*PMiNi*Pl  * 

FMLT*ABS(OUi*M/UI>PK  ) 

(■OMHYsni'MN^p^Ll  *  (  OU^Mt.'4-l)nyN  ) 

AtPHA»  (PNP/(nr>PMY«lll)t  K  )*PPLS/'JIJPK>/K^t(» 

Pfc  T  A  ••PPL  S*  (  ALPHA*  t  ,tf /MOPK  ) 

PNGZ«»-MP/r<PLS 
nNOM«AI.PHA*RWP*Ht  TA 

PCP^NP  1  f  P\P/1/NOM*U(iPk  ) 

rPMM  H ■  A (  (>G(»<C»MN) 

f  GZL»(HE  T  A*f  I  .i./hcnmiy.j  #n)  /  HNllM)/f  { t  NU/*Cr*Nl  h*»np**  (PNGZ-t  ,P)  )* 
t  (  WNP*OI)Pk  *r  »:om)  ) 

BGZ*fc  APlHbZL 1 

C 

to  »HW«PPk-R 

r.P«i  } 

Tt  (W.GT.k/j  &P«rRPc-KWPPLS*GR*{R-«Z)/PPlS 

PpaHf<K/(A».PMA*PHWAHt  TA)  *U0PK 

nOP»GP#H» 

PPKbWP**  \  tf  e, 

PbPai  ,F5 

1*  CT.CiT.J.vi)  («n  Til  )3 
C 

11  t  F  (T.t*).TfPU)l  (ill  re  l  y 

c 

a«-i ,2e-3 

C  •  a  LUG  f  •  A  /  f  (ll’PI*  •  A  )  )/(K7(j.RpK) 
pa  tni)Pn-A  )  AtxP(-C*wPK  ) 

C 

12  wFl TaA*H#E*PfC*W) 

JF  n.LE.P.2)  OoPaaFLT 

TF  (T.GT.O.?)  0D*tatHPLT*(l.-T)A0iJB*CT«O,2)  )*1  ,25 
C 

13  TTOLt)«T 
PtTUKN 

c 

14  nKUHM^n 

II ad*  F0R*A1(»  rPPO»|  TIMt  PU8T  PE  GREaTE«  THAN  ZFROM 
8T0P 
FNO 


_ US'  . 
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